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Abstract. In this paper we present a probabilistic algorithm to compute the coefficients of modular forms 
of level one. Focus on the Ramanujan's tau function, we give out the explicit complexity of the algorithm. 
From a practical viewpoint, the algorithm is particularly well suited for implementations. 



1. Introduction and Main Results 

In the book[10], Couveignes, Edixhoven and others describe an algorithm for computing coef- 
ficients of modular forms for the group SL2(Z), and Bruin [6] generalizes the method to modular 
forms for the groups of the form Ti(n). Their method leads to a polynomial time algorithm for 
computing coefficients of modular forms. But, efficient ways to implement the algorithm and ex- 
plicit complexity analysis are still being studied. Working with complex number, Bosman's explicit 
computations shows the power of this new method. As one of the applications, he largely improved 
the known result on Lehmer's nonvanishing conjecture for Ramanujan's tau function. Following 
Couveignes's idea[8], we give out a probabilistic algorithm, which seems to be more suitable to do 
complexity analysis. Instead of using Brill-Noether's algorithm, we work with the function field of 
the modular curves, using Hefi's algorithm to compute the jacobian of modular curves. 

The discriminant modular form is denoted by, 

CO CO 
71=1 71=1 

where z € H, q = e 2nlz . 

Let I be a prime, the modulo £ Galois representation associated to A(q) lies in the Jacobian 
variety of the modular curve X\ {£) , which can be written as a finite intersection 

V t = P| kerCZl-T^Ji^XQF]) 

l<i<(£ 2 -l)/6 

where Tj, i = 1, . . . , are Hecke operators. Indeed Ve is a group scheme over Q of order £ 2 , which 
is called Ramanujan subspacef ]. Once Vg mod p have been computed for sufficient many small 
primes p, then V^/Q can be recovered. Given Vg, the calculation of t(P) mod £ for large prime P 
is easy. 

Though out the paper, p , £ denote different primes and g = - 7 ^ is the genus of X±(£), 
unless otherwise stated. The complexity of single group operation in the Jacobian will be deonted 
by 0(<7 W ), where a; is a constant in [2,4], see[16][15], and e is any constant bigger than 0. The main 
results of the paper are as follows, 

Theorem 1.1. Let £,p > 7 be different primes. Given the Zeta function of modular curve Xi(£)f p 
then, Vi mod p can be computed in time 0((£ 4+2uJ log 2 p) 1+e ) . 
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Corollary 1.2. Let 6 be a constant such that for any prime £ > 7, there exists nonconstant function 
h : Vg — > Pq, satisfying the height of coefficients of polynomial P(X) := rTaeV^-O^ ~~ h(a)) is 
bounded above by 0(£ b ), then 

(1) Ve can be computed in time 0(£ 4+2w+2b+e ). 

(2) For prime p, r(p) can be computed in 0(log 4+2w+26+e p). 

Remark 1.3. Using Khuri-Makdisi's algorithm, the constant uo can be 2.376. Using Hefi's algorithm, 
u is known to be in [2,4] only. The constant b is less than 16, showed in[10], which seems to be 
rough. Our computation suggests that, b is something like 3. 

Notation: The running time will always be measured in bit operation. Using FFT, multiplica- 
tion of two ra-bit length integers can be done in 0(n 1+e ) time. Multiplication in finite field ¥ q can 
be done in 0(log 1+e q). 

The paper is organized as follows. Section two provides some necessary background on computing 
the plane model and the function field of modular curves, computing isogeny of elliptic curve over 
finite field are also recalled. A better bound on the generators of the maximal ideals of Hecke 
algebra is proved, which is used to reduce the complexity of the algorithm. 

Section three contains the application of HeB's algorithm to the computation of the jacobian of 
modular curves over finite fields. Where, we give out methods to find the correspondence between 
the places of function field and the cusps of modular curve, and compute the action of Hecke 
operators on places of function field. 

The main algorithm is given in section four, including complexity analysis. 

The final section is concerned with some real computations of Ramanujan's tau function. All of 
our computation are based on Magma computational algebra system[3]. 

2. Function field of modular curves 

In this section we study the plane model and the function field of modular curve X\{£). Let T\{£) 
be the congruence subgroup of SL2(Z), % be the upper half complex plane, the so-called modular 
curve over C is defined as 

Y 1 (£) = r 1 (£)\n 

One can add cusps P 1 (Q) to Y\{£) to compactify it and obtain the modular curve 

x 1 (£) = r 1 (£)\n* 

this complex algebraic curve is defined over Q, denoted by Xi(£)q, moreover, for £ > 5, X\(£) 
has natural model over Z[l/£|. Let K be a number field, then K- valued points of Y\{£)q can be 
interpreted as following, 

Y^MK) = {(E,P) : E/K,Pe E[£](K) - {0}}/~ 

where E is an elliptic curve over K, P is a ET-rational point of order £, and (Ei,P\) ~ (^j-fb) 
means that, there exists a iT-isomorphism (f> : E\ — > E2, such that 4>(Pi) = Pi- Such a mod- 
uli interpretation implies a way to obtain the plane model of modular curves, as the following 
proposition[2], 

Proposition 2.1. Suppose that £ > 4. Then every fT-isomorphism class of pairs (E,P) with E 
an elliptic curve over K and P G E(K) a torsion point of order £ contains a unique model of the 
Tate normal form 

E{b,c) ■ V 2 + (1 " c)xy -by = x 3 - bx 2 , P = (0, 0) (2.1) 

with c S K,be K*. 

Thus, points of X\{£) can be represented as pairs (b,c) in a unique way. The £-th division 
polynomial gives out a polynomial in b and c, which defines a plane curve birationally equivalent to 
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X±(£), after carefully change of variables, the polynomial of the plane curve becomes much simpler. 
We use the table listed in Sutherland's homepage, for example the plane model of Xi(19) is, 

f(x, y) =y 5 - (x 2 + 2)y 4 - (2x 3 + 2x 2 + 2x- l)y 3 + (x 5 + 3x 4 + 7x 3 + 6x 2 + 2x)y 2 

- (x 5 + 2x 4 + 4x 3 + 3x 2 )y + x 3 + x 2 

where 

x(x + y)(y — 1) x(y — 1) 

r = l + - ; y ' yy — -, s = l + - ^ 1 2.2) 

(x + l)(x z — xy + 2x — y z + y) (x + l)(x — y + 1) 

and 

c = s(r — 1), b = rc (2-3) 
so the function field of X\(19) over Q is, 

QCXi(19)) = Q(6,c) = Q(x)[y]/(/(x,y)) 
For simplicity, throughout the paper, the affine model of modular curve X±(£) will be denoted as 
/(x, y), and r = , s = g^'yj , c = s(r — 1), b = rc are the coordinate transformations. Notice 

that, pair (b,c) satisfying, -E(b, c ) is an elliptic curve with point (0,0) of order £. Cusps of X\{£) 
corresponding to those pairs (b, c) such that the j-invariant j(E^ b c ^) = oo. For £ an odd prime, 
X\[€) has £ — 1 cusps, half of which are in Xi(£)(Q), the other cusps are defined over the maximal 
real subfield of Q(C^)- Interpreted in the function field, every Q-rational cusp corresponding to 
a degree one place of Q(Xi(£)), denoted as 0\, . . . ,0(^_i)/2, the Q(<^)-cusps corresponding to a 
degree place of Q(Xi(£)). It's easy to get these places after writing down the exact expression 
of j(-^(b, c )) m variables x,y. As the above example, places corresponding to Q-rational cusps are, 



Oi 

o 2 
o 3 
o 4 
o 5 
o 6 
o 7 

Os 

o 9 



y 4 y 3 (-x 3 - x 2 + x - 1) y 2 (-x 3 - 2x 2 - 2x - 2) 2y 2x - 1 

x 4 + x A x 4 + x^ x^ + ar x x 

y 3 y 2 (x + l) 
1 hy), 

X X 

y 4 y 3 (x-l) y 2 (-x 3 - 3x 2 - 3x - 2) 2y x-1 

x 4 + x J x 4 + x J x- 3 + x z x x 

, , 2y 3 2y 2 

» + Hy), 

X X 

y 4 y 3 (-2x 3 - 2x 2 - x + 1) y 2 (-x 3 + 3x 2 + 5x + 2) y(x - 2) x + 1 . 

X — T" 1 j T q "I - 7 o H I - o o ~~ f~ H~~ 

X* + X° X^ + X J X J + X' X X 

1 y 4 y 3 (x 2 + l) y 2 (2x 2 + 3x + 2) y(-x 2 - 3x - 5) 

-,-^4+1 , 

X X 4 

1 y 4 y 3 (-x 2 -2) i y 2 (-2x 2 - 2x - 2) i y(x 2 + 3x + 7) i 1 
x x a 

i y y y 



+ — g L + — s - + - 5 ~ + 

X° X 4 X z X 



rj* ' CC*^ 



+ 1), 



where places are represented as prime ideals of orders of the function field Q(Xi(£)). 

It's well known that[9], modular curve X\{£) has good reduction at those primes p \ £, the 
reduction curve is denoted by Xi(£)^ . Having the nonsingular affine model of Xi(£), it's easy to 
calculate the affine model for Xi(£)^ p , hence the function field of X\{£)^ . For simplicity, the plane 
model (resp. coordinate changes and points ) of X\ (£)v p , will still be denoted by /(x,y) (resp. r, s 
and (b, c)), and Oi,i & {1, . . . , ^-} denote the F p -rational cusps of Xi(£)^ p , which are the reduction 
of Q-rational cusps of X\ (£) . 
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The Ramanujan subspace Vi mod p is a subgroup scheme of Jacobian variety (denoted by J\ [£)w p 
) of Xi(£)f p , defined as following 

where Tj, 1 < i < are Hecke operators, the number follows from[19]. In fact, the Hecke 

algebra T = Z[T n ,n G Z + ] C End(Ji(£)) is a free Z-module of rank g = — 5 ^ 7 ^ , after representing 
each Hecke operator as a matrix, the generators can be extracted from T\, . . . , Tjv2_i)/6 by solving 
linear equations. For example, when the level £ = 17, the Hecke algebra T is equal to ZTi + . . .+ZT48 
as a Z- module, which can be changed to ZTi + ZT2 + +ZT3 + ZT4 + ZTg as a free Z-module, so there 
are fewer Hecke operators and isogeny of lower degree needed to take account. But, in practice, 
we can do much better, notice that the goal is to find nonzero element in J\ (£)■$[£], such that 
after the action of Tk — r(k),k > 1 it will be zero. Assume there is an element D G Ji(^)f p M 
satisfying (T 2 — t(2))(D) = 0. If one has the relations Tk — r(fe) = <pk ■ (T2 — t(2)), for some 
endomorphism tfik G End(Ji(£)p p ) in advance, then (Tk — r(k))(D) is equal to zero automatically 
(so D G Vi). The action of Tf, — r(k) on Ji(£)w p [£] can be represented by a matrix over finite 
field ¥1, the existence of 4>k is equivalent to the existence of some matrix Mjt over such that 
T k - r(fc) = M k ■ (T 2 - t(2)). For example, when the level £ G {13, 17, 19, 23, 29, 37, 41, 43}, for any 
k > 3, Tk — r(fc) = Mfc • (T2 — t(2)), for some matrix Mk over finite field ¥1, while I = 31, T3 — r(3) 
has this property. In fact, we proved the following proposition, 

Proposition 2.2. Let I be a prime bigger than 3, 5 2 (ri(£)) is the space of cusp modular forms 
of weight 2 level I over C. Let T = Z[T n \ n > 1] C End(5 2 (ri(^))) be the Hecke algebra and m 
is the maximal ideal generated by £ and the T n — r(n) with n > 1. Then m can be generated by I 
and the T n - r(n) with 1 < n < [^±i] . 

Proof. Let ^ := [2^1] and R := T ® z F f , then is an Artin ring, which can be decomposed as 

R = II R e 

p 

where p run through all maximal ideals of R, and R p is the localization of R at p. 

Let m be the image of m in R. Then, it is enough to prove m can be generated by the T n — r(n) 
with n < hi. Or equivalently, to show, any Tk — i~(k), k > bi can be represented as 

be 

T k -r(k) = Y,A-(T l -r(i)) (2.4) 
i=i 

where Ai are operators in R. 

Let 52(ri(^);F^) := 52(ri(^);Z) <8>F^, which is a i?-module, decomposed as 

S 2 (T 1 (£);¥ e ) = HS 2 (T 1 (£);¥ e ) p 

p 

where p run through all maximal ideals of R, and S2(T\{£);¥t)p is the localization of 5 2 (ri(£);F^) 
at p, which is a i? p -module. To show 2.4, it's enough to show for each p, the action Tk — r(k), k > be 
on S2(Ti(£);¥i) p can be represented as 

bi 

T k -r(k) = Y,Bi-(T i -T(i)) (2.5) 
i=i 

where Bi are operators in R p . 
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The maximal ideal p of R corresponds to a Gal(F^/F^)-conjugacy class [/] of normalized eigen- 
forms in S2(Ti(£);¥i), they are newforms of level I. Thus the localization S2(Ti(£); ¥e)p is a vector 
space spanned by these newforms, so R p is a field isomorphism to the field F^(/), which is generated 
by the coefficients of /. 

Let a,i(f) denote the i-th coefficient of /, to show 2.5 it is enough to show a,k(f) — r(k) is always 
equal to zero, or there exists at least one nonzero element in the aj(/) — r(i) with 1 < i < bg. 

If / is congruence to A (5) mod £, then afc(/) — r(k) is always equal to zero. Now suppose / is 
not congruence to A(q) modulo £ and a,i(f) — t{i) = for 1 < i < bg. The Proposition^ 10(b) of [11] 
together with Theorem3.5(a) of [1] implies that / comes from a level one newform of weight k\ with 
fci < 2£. Since r{p 2 ) = r(p) 2 - p u and a 2 (f) = a p (f) 2 - p kl ~ l , we have p kl ~ 12 = 1 mod £ for all 
primes p < bg_. Let p the minimal primitive root modulo £ (which is in 0(log 6 p) if assume GRH), 
so fci - 12 = mod £-1. Let A(q) = 1 be the Hass invariant, then fx := A{q) kl ~ 12 ^ l ~^ A(g) is 
a newform of level one weight fci, and dk(fi) = dfc(/) for all k < b(_. Since bp > Yt 1 the theorem 
of Sturm[19](or see [7]) tells that f\ is congruence to /, that is A(q) is congruence to / modulo £, 
which is a contradiction. □ 

Since £ and a subset of {T2 — r(2), . . . ,Td — r(d)} may also suffice to generate m, define the 
optimal subset as follows, 

Definition 2.3. Let S be a set of positive integers, such that £ and T n — r(n),n G S generate m. 
The set S is called optimal if m&x{p \ p\n, n £ S} is minimal among all the 5. 

In practice only those operators in the optimal subset need to consider, which will faster the 
algorithm. 

Let n > 2 be a prime not equals to £, and Q a point of Xi(£)f p (¥ p ) represented by (E( b ^, (0, 0)), 
the computation of T n (Q) comes down to the computing isogeny of elliptic curves over finite field, 

T n (E {biC) ,(0, 0)) = ^(E/C, (0, 0) + C) 

c 

where C runs over all the order n subgroups of En, jC \. The following result about the complexity 
of computing n-isogeny is in[5], when n is small compare to the characteristic of the field, 

Proposition 2.4. Let ¥ q be a finite field of characteristic p, n a prime not equals to p, E and E 
are elliptic curves over ¥ q , such that there exists an isogeny (p : E — >■ E of degree n, then (j> can be 
computed in 0(n 1+e ) multiplications in the field F„. 

There are two things should be mentioned, firstly, the original elliptic curve £Vb lC ) is defined 
over some finite field F, but the isogeny may lies in extension field of F, the exact definition field 
of isogeny can be obtained by computing the n-th classical modular polynomial $> n (X,Y), and 
solving the equation $ n (X, j(E(^ c \)) = 0, the extension degree is less or equals to the degree of 

$ n (X, j^Ei^a)). Secondly, in general the isogenous curve E is not in Tate normal form 2.1, using 
the map (f>, or rather the point c/>((0, 0)) E E, E can be transformed to Tate Normal form after 
some coordinate changes, this gives a new point of modular curve Xi(£)y p - 

3. Computing the jacobian of modular curves 

Let ¥ q be a finite extension of finite field F p , Xi(£)p q is the base change of X\(£)f p to ¥ q . One 
of the important task is computing the Jacobian Ji(£)w q of modular curve Xi(£)f q . For general 
curves, we already have polynomial time algorithms to perform operation (addition and subtraction) 
in their Jacobians[14][20][15][16]. For Jacobian of modular curve, Couveignes uses Brill-Noether 
algorithm to do the computation [8], while Bruin uses Khuri-Makdisi's algorithm[G], in this paper, 
we choose Hen's algorithm, with the advantage that it's easy know the relations between points of 
modular curves and place of its function field. 



(i 



JINXIANG ZENG 



Let's recall the main idea of HeB's algorithm, for the detail see[15]. Let K = ¥ q (x)[y]/(f(x,y)) 
denotes the function field of modular curve X±(£)w „• There are isomorphisms, 

Ji(£)(¥ q ) * Pic°(Xi(£) F ,) = C1°(K) 

Notice that, there is some calculation behind the second isomorphism, i.e. computing the change 
of representations. HeB's algorithm is based on the arithmetic of function field K. Let P be the set 
of all places of IK, and S be the set of places satisfying vp{x) < 0. The ring of elements of K being 
integral at all places of S resp. P — S is denoted by Os resp. O s , they are Dedekind domains, 
which are called infinite and finite order of K. Their divisor groups are denoted by Div(0,g) resp. 
Div(0 5 ). Place of K corresponding to prime ideal of O s or O s ', in fact Div(K) can be decomposed 
as [15] 

Div(K) ^> Div(0 5 ) x Div(C 5 ) 

As O is a Dedekind domain, its prime ideal can be represented by two elements of 1K X , which can 
be normalized as f\{x) = x d + a^-i 3 ^ -1 + . . . + a\x + 1, and f2(x, y) = y m + b m -\{x)y m ~ 1 + . . . + 
bi(x)y + bo(x), where ai,0 < i < d — 1 belong to the constant field ¥ q , and bi(x), < i < m — 1 are 
elements of F g [x], with degree less than d. So the point set corresponding to a prime ideal p of O s 
can be computed as follows, let fi(x), f2(x, y) be the normalized two generators of p, compute the 
roots of fi(x) = 0, denoted by Xi, 1 < i < d, for each root Xi, compute roots of f2(xi,y) = 0, denoted 
by yij, 1 < j < m, then the point set corresponding to p is {(xj, yij) | 1 < i < d,l < j < m}, notice 
that every point (xi,yij) satisfying the plane model f(xi,yij) = 0. Conversely, given a point set 
{(x{, y^) | 1 < i < d, 1 < j < m}, the corresponding prime ideal of O s can be computed as follows, 
the first generator is clear, f\{x) = \\/ l= i{x — x^), the second one need some calculation, assume 
bi(x) = Xlfe=o c ik xk 1 < i < m — 1, where are parameters belong to ¥ q , after interpolating 
the points (xi,yij), it's a linear equation of md variables, solving this equation gives the second 
generator. 

So, a degree d place P of K = ¥ q (x)[y]/(f(x, y)), corresponding to point set, denoted by 
{(xi,yi),l < i < d}, which forms a completely Gal(F d/F g )-conjugate set. Using the coordinate 
transform formulae (2.2), for each point (xi,yi), the corresponding point of X\(l)^ v in the form 
(£/(6 ijCi ), (0,0)) is clear, where En, iiCi ) is an elliptic curve in Tate normal form, and (0,0) is a point 
of order I. As discussed in section 2, the action of Hecke operator T n on each point (E 1 ^.^), (0,0)) 
leads to a sequence of elliptic curves. Using the inverse transform formulae these curves give the 
point sets in the affine curve f(x,y) = 0. Further, one has places of the function field K, 

Definition 3.1. Let P be a place of K = ¥ q (Xi(£)), with degree d , T n is a Hecke operator. T n (P) 
is defined to be the divisor of K, corresponding to the point set ^^ =1 T n (E( fe . c .), (0,0)), which is 
effective of degree ip(n)d. If D is a divisor of IK, decomposed as X12=i a *-^' where Pi are places, 
then T n (D) is defined as X^Li a iT n {Pi)- 

Let Dq be a fixed degree one place of IK, which will be served as an origin, then every element of 
Cl° (IK) is represented by D — gDo , where g is the genus of IK, and D is an effective divisor of degree 
g. Addition in C1°(1K) means that, given A, B effective divisors of degree g, find an effective divisor 
D of degree g, such that D — gDo is linearly equivalent to A — gDo + B — gDo. The complexity of 
this operation can be found in [15], 

Proposition 3.2. Notations as above. There exists a constant uj € [2, 4] such that, the divisor D 
can be computed in 0(g UJ ) multiplications in the field ¥ q . 

Remark 3.3. We have not yet seen the precise value of uj. However, using Khuri-Makdisi's algorithm 
the complexity of single group operation is known i.e. uj = 2.376, when fast algorithms for the linear 
algebra is used. 
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3.1. Searching an ^-torsion point. One of the main steps in computing the Ramanujan space 
Ve mod p is to find an ^-torsion point in Ji(£)(¥ q ), in order to do that one has to work with some 
large enough extension field of ¥ p , such that Ji(£)(¥ q ) contains ^-torsion point. A direct way to 
get such a point is, pick a random point Qo in Ji(£)(¥ q ) and then compute Q\ := NgQo, where 
£ e Ng = #Ji(£)(F g ) and I \ Ng. If Q\ is nonzero and ^-torsion, then it's good, if it's not, try 
Q2 := £Qi, check again, after several steps, we get a nonzero ^-torsion point. As #Ji(£)(F q ) is 
bounded above by q 9 , using fast exponentiation, the running time of getting an ^-torsion point 
is about log(g 9 ) • O(^) = 0(g l+UJ logq) multiplications in the field ¥ q or 0(g l+w log 2+£ q) bit 
operations. The computation is costly, due the huge factor Ng, which is nearly q 9 . However, we 
can faster the calculation by introducing some tricks. 

For each newform / in S2(Ti(£)), there is an associated abelian variety Af, and ri/Af is an 
isogeny decomposition of Ji(£), where / runs through a set of representatives for the Galois conju- 
gacy classes of newforms in S2(Ti(£)). Let fi be the newform, which is congruence to A(q). Then 
the Ramanujan space Vg land inside Af e . We can detect Af v since we know the characteristic poly- 
nomial of Frobenius acting on Vi mod p is X 2 — r(p)X +p n mod £, which should be a factor (over 
F^pT]) of the characteristic polynomial of Frobenius acting on Af e mod p. Now let A' := Af be 
the product of abelian varieties, except Af e . The minimal polynomial of Hecke operator T 2 acting 
on A' is denoted by P 2 [X] € Z[X\. Now for any random point Q in Ji(£)(¥ q ), Qi := P 2 (T 2 )(Q ) 
will be a point in Af e (¥ q ). Similarly, let Ng be the prime-to-£ part of #Af e (¥ q ), by computing 
NiQi, we can obtain an ^-torsion point of Ji(£)(¥ q ). Since the calculation of P2(T 2 )(Qo) is easier 
and Nn is smaller, which is bound above by q to the dimension of Af e , we can faster the algorithm, 
especially when the dimension of Af e is small compare to g. 

We now explain how to compute #Ji(£)(¥ q ), 

Lemma 3.4. (Manin, Shokurov, Merel, Cremona). For £ a prime and p another prime, the Zeta 
function of X\(£)-g can be computed in deterministic polynomial time in £ and p. 

Given this Zeta function one can easily compute #Ji(£)(W q ), for q is a power of p. 
The following definition[15] is very useful in our algorithm, 

Definition 3.5. Let A be a fixed place of degree 1, and D a degree zero divisor, then Do is called 
the reduction of D along A, if Dq = D + dA with d > is minimal. 

It's easy to see the reduction divisor D$ is effective and unique. 

3.2. Distinguishing the rational cusps. The Hecke operator T n ,n € Z + is defined over Q, so 
it takes Q-rational points to Q-rational points, the action of T n on Q-rational cusps is needed in 
the algorithm, it seems that, there is no an easy way to know the action directly from the place 
Oi,i € only. But, a coset representatives of Q-rational cusps of X±(£) is clear, i.e. 
{ T> 2' • • • ' (1-1)72 }' anc * f° r P r i m e n / £, T n (—) = — + n=, where nfn is the class of integer nm in 
(Z/£Z)*/{±1}. So knowing the 1-1 correspondence between {1, 1 < % < ^} and {O h l<i< ^-} 
, leads to knowing the action T n (Oi). In general, it's not easy to know the correspondence[17]. Our 
strategy is to reduce the problem to finite field, as follows, choose a prime p, such that ffJi(£)(¥ p ) 
has a small factor d, let g be the genus of the function field ¥ p (Xi(£)). Fix a cusp Oj served as 
origin, as discussed above, let D — gOi be a degree zero random divisor of order d, in general 
D doesn't contain cusps, if it does, try a new one. Now, assume Oi corresponds to ^ for some 
m 6 {1, . . . , ^p} and Oj corresponds to == for some j € {1, . . . , —^}- Then, for prime n, compute 
D n := T n (D — gO{) by the assumption as follows, 

D n = T n {D) - gT n (Oi) = T n (D) - 5 (0, + nOj) 

if D n is not of order d, then the assumption is wrong, replace Oj by another — or = by another 
Oj, try again. The complete correspondence will be detected after several tries. In fact, the 
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correspondence is known up to cyclic permutation, but it will be enough for our algorithm. There 
is one more thing should be noticed, the degree divisors of the form Oj — Oj, 1 < i,j < ^s- 
generate a subgroup of Ji(£)(¥ p ), which is called cuspidal subgroup, so the chosen factor d should 
not be a multiple of the order of cuspidal subgroup. 

As the example in section 2, the correspondence between {0\, . . . , 0$} and {|, . . . , ^} is 

111111111 

156743289 

Ox 2 3 4 O s O e 7 8 9 

After the above preparation, we can now estimate the complexity of computing T n (Q), where 
Q is a degree zero F 9 -divisor of the form YlT=i a *^« ~~ ff^> ^ * s fi xe d to be the cusp 0\ and n is a 
prime. 

As mentioned above, first compute the point set of the degree d{ place Pi, pick one of them forms 
a point (i£(6 i)Ci ), (0,0)) on the modular curve, which is defined over ¥ d t . Denote the factorization 

of Q n (X, j{E(p^ Ci })) over F^JX] as Y\i=iFj(X), with Fj(X) irreducible of degree fj. So for each 
root of Fj(X) = 0, there is an isogeny of degree n defined over F d t f 3 , computing this isogeny takes 

0(n 1+e • (lo g( f^) 1+£ ) = Odndifj log q) 1+e ) 

bit operations. 

Notice that, isogenous curves corresponding to the roots of Fj(X) forms a Ga^F^/j /F^)- 
conjugate set, so compute any one of them will be suffice. So the complexity of computing 
T n ((E ktCi ,(0,0))) is about 

h 

Q:=^0((nd i / i log( 7 ) 1 + e ) 
i=i 

Since Ylj=i fj = nJl ~ 1> ^ ^ s bounded above by 0((ndin log <?) 1+e ). 

The complexity of computing T n (Q) is C := X)2=i Since Y^Li = 9i @ ls bounded above 
by 0{(n 2 g log q) 1Jre ) bit operations. 

Notice that after the action of T n , the divisor T n (Q) becomes complicated, so reduce it along 
O before going into the next Hecke operation. The complexity of reduction is about 0(logn • 

log 1+e q) bit operations, where the factor log n means the times of reductions need to perform. 

4. Computing the coefficients of modular forms 

After some modification, the simplified algorithm proposed by in[8], tends to be very efficient to 
compute the Ramanujan subspace Ve mod p. We can even give out a complexity analysis of the 
algorithm. The main idea of the algorithm is, pick some random points in Ji(£)(W q ), make them 
^-torsion points, then project them to the space Ve mod p by the Hecke operators. Using Chinese 
Reminder Theorem, V^/Q is recovered from sufficient many Ve mod p. 

4.1. Computing Ramanujan subspace modulo p. The characteristic polynomial of Frobenius 
F p acting on Ramanujan space Ve mod p is Fi 2 — r(p)F p + p 11 mod £, which gives the definition 
field of elements of Ve mod p, denoted by F g . It will be nice to have sufficient many primes p such 
that F g is small, in fact, the following holds, 

Lemma 4.1. Let A := r(p) 2 — Ap 11 , then we have the following results, 

(1) If A = mod £, then (F p ) r is identity, for some positive integer r\£{£ — 1). 

(2) If A mod £ is a nonzero square, then (F p ) r is identity, for some positive integer r\(£ — 1). 

(3) If A mod £ is not square, then {F p ) r is identity, for some positive integer r\(£ 2 — 1). 
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Proof. The Frobenius endomorphism F p acts on Vi mod p via a 2 x 2-matrix with characteristic 
polynomial F 2 — r(p)F p + p 11 = 0. If the matrix has a double eigenvalue, let r be order of the 
eigenvalue in Wi, then the order F p is r resp. £r, while F p is diagonalizable resp. not diagonalizable. 
If the matrix has two different eigenvalues in F^, then the order of F p is the least common multiple 
of the orders of these two eigenvalues, which is a factor of £ — 1. If the matrix has two conjugate 
eigenvalues in ¥ p — F^, then the order of F p equals to the order of the eigenvalue, which is a factor 
of £ 2 - 1. □ 

So, for at least half of the primes p, the definition field of Vi mod p is an extension field of ¥ p of 
degree in 0(£), denoted by ¥ q . Set the base field to be ¥ q , then the space Vi mod p is a subgroup 
of Ji(£)(¥ q )[£]. Notice that, elements of Ji(£)(¥ q )[£] may lie in extension fields of ¥ q . 

For every integer n > 2, the characteristic polynomial of T n acting on S2(T\(£)) is a degree g 
monic polynomial belongs to 1>[X], denote it by A n (X), which can be factored as 

A n {X) = B n (X)(X - T (n)) e " mod£ 

with B n (X) monic and B n (r(n)) F^, the exponent e n is > 1 due to the theorem of congruence 

of modular forms (Thm2.5.7 ofjlO]). We call ir n : J 1 (£)(¥ g )[£] — > Ji(£)(¥ q )[£] the projection map, 
which maps an ^-torsion point D £ Ji(£)(¥ q )[£] to an ^-torsion point B n (T n )(D) of Ji(£)(¥ q )[£], 
and maps bijectively Vi mod p onto itself. If E := ir n (D) ^ 0, the minimal exponent d n is defined 
as the nonnegative integer satisfying, 

(T n - r(n)) d -(E) + and (T n - r(n)) d - +1 (E) = 

since (T n — r(n)) en (E) = 0, then < d n < e n . Let 7r n be the composite map of ir n and (T n — T(n)) dn . 
Let S be a optimal set as defined in Definition2.3, then for a random element D € Ji(£)(¥ q )[£], 
tts(D) := W n& g^n{D) is an element in Vi mod p. Knowing the characteristic polynomial of 
Frobenius F p acting on Vi mod p will help us to find another element of Vi mod p with lower 
cost in some cases. Namely, if the characteristic polynomial has two different eigenvalues then 
for any nonzero element E € Vi mod p, F p (E) and E is linearly independent, so Vi mod p = 
¥iE + ¥i(F p (E)), if the characteristic polynomial has same eigenvalue, then F p (E) = \E, where 
A G F| is the eigenvalue. In this case, choose a new random point D\ of Ji(£)(¥ q )[£], compute 
irs(Di), then check whether tts(D) and tts(Di) is linearly independent. If it's true then success, if 
it is not, try another random point Di- For two nonzero random points D±,D2 € Ji(£)(¥ q )[£], the 
possibility of tts(D\) and ^5(1)2) is linearly independent is nearly j, so after several tries we will 
get a base of mod p. 

4.2. Complexity analysis. Now we can prove Theoreml.l. Let p be a good prime, first compute 
the definition field of Vi mod p, which is an extension of F p of degree in 0(£), denoted by ¥ q , the 
complexity of getting an ^-torsion point of Ji(£)(¥ q ) is about 0(g 1+UJ log 2+e q) bit operations. 

The second step is to compute the action of Hecke operators. Let B n (T n ) = T d + ad~iT d ~ l + 
. . . + a\T n + ao, where ai G F;, the computation of B n (T n ) acting on a divisor Qo := D — gO can be 
done as follows, firstly, compute Qi+i := T n (Qi), then reduce it, then compute Qi+2 '■= T n (Qi + i), 
and so on. Given Qi, the complexity of computing and simplifying Qi+i := T n (Qi) is about 

0{(n 2 g log + 0(^ log n log 1+e q) 

As d is in 0(g), computing and simplifying Qo, ■ ■ ■ ,Qd takes 

d(0{(n 2 glogq) 1+£ ) + 0( 5 W lognlog 1+e g)) = 0{(n 2 g 2 log q) 1+e ) + 0(g 1+u> lognlog 1+e g) 

Given Qi, as <2j € F^ the computation of aj • Qi takes 0(log^ • g w log 1+e q), so the computation of 
Q d + a d ~iQd-i + • • • + aiQi + a Q takes 

Oidg" log£log 1+e q) + 0(log d • g" log 1+e q)) 

As d € 0(g), it's bounded above by 0(g 1+UJ log£log 1+e q). 
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In summary, the complexity of B n (T n )(Qo) will be 

0({n 2 g 2 log q) 1+e ) + 0(g 1+u> log n log 1+e q) + 0( 5 1+tJ log £ log 1+e q) 

which is C n := 0(((ng) 2+€ + g 1+UJ (logn + log^)) log 1+e q) . The total complexity of computing 
ns(Qo) = Hues B n(T n )(Q ) is C := J2nes C n- As n G <S,n < £ and |<S| < £, C is in 

0(((^V) 1+e + g 1+LU £ 1+e ) log 1+e g) 

since g = £f± and g € 0(/), C is in 0({£ s+e + £ 4+2w+e ) log 1+e p). As w is in [2,4], so C is 
O (£ 4+2a)+€ log +e p) . The complexity of getting an ^-torsion point is in 0(£ A+2uj log 2+{E p). The com- 
plexity of checking linearly independence of two elements in Vg mod p is about 0(^log 1+£ g) = 
0(£ 2+UJ+e log 1+e p), so the complexity of computing a base for Vg mod p is in 0(£ i+2ui+e \og 2+t p). 

Remark 4.2. (l)As showed in above, the algorithm takes the main effort to get an ^-torsion point. 

(2)If without Proposition. 2, the complexity will comes up to 0((£ w+e + £ 5+2ul+e ) log 1+€ p) + 
0(£ 4 + 2a, + e log +e p), the dominated term is computing the action ns- 

4.3. Computing Ramanujan subspace. Fix a Q-rational cusp 0\, denoted as O, which will be 
served as the origin of the Jacobi map. For every element x € Vg, let D be the reduction of x along 
O, i.e. D = x + dO, where D is an effective divisor of degree d, denoted by D = Q\ + . . . + Qd and d 
is called the stability of x, denoted by 6(x). Choose a rational function h(x) £ Q(Xi (•£)), which has 
no pole except at O, and define a function G : Ve — > Q as G(x) = h(Qi) + . . ■+h(Qd) € Q. From the 
uniqueness of D, one has G(a(x)) = a(G(x)) for any a € Gal(Q/Q), which is very important. Let p 
be a prime ideal of Q(x) over p, we hope the uniqueness property still holds in the reduction world, 
which means the following: let x (resp. O, D, Qi) be the reduction of x modulo p (resp. O, D, Qi). 
The reduction of x along O is denoted by E = x + 0(x)O, where 9{x) < 6{x). If 9(x) = 9(x), then 
by the uniqueness property E = D ( =Q± + . . - + Qd )■ As we don't know the value 9{x) in advance, 
a criterion to determine when 9{x) = 9(x), Vx G V; is needed. In theoretically, Bruin[6] proved 
that, for at least half of the primes smaller than £°^\ the following holds: 9(x) = 9(x), Vx G Vg, 
we call them good primes. In practice, since 9{x) is less or equal to the genus of Xi(£), so if 9{x) 
is equal to the genus for every x £ Vg mod p, then 9(x) = 9{x) automatically for all x G Vg . 
Now the polynomial 

P(X):= H (X-G(x)) 
xeVt-o 

is a polynomial in Q[X] of degree £ 2 — 1, whose reduction modulo p is exactly the polynomial 

P(X) := J] (X - G{x)) 

xSVf mod p—O 

where G is the reduction map of G. 

Given Vg mod p = ¥ge\ + the computation of P(X) comes down to perform £ 2 additions 
in the Jacobian, with a complexity 0(£ 2 £ 1+2u log 1+e p) = 0(£ 3+2u ' log 1+e p). 

Since the height of P(X) is expected to be in 0(£ b ), for some integer b < 16 ([] ^Theorem 11.7.6). 
Using the fact, there exists a constant c such that 

p > cexp(L) 

p<L, prime 

to recover P(X) from P(X)'s, it will be suffice to take the upper bound L of good primes to be 
0(£ b ), so the complexity of computing P(X) will be 

0((£ 4+2 " log 2 p) 1+t ) = 0{{£ A+2uJ L 2 ) 1+ *) = (£ 4+2 - +2fe+e ) 

p<L, prime 
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In practice, choose the function h{x) € Q(-Xi (•£)), such that the degree is equal to the gonality of 
the curve X\(l). For t < 40, h{x) have been computed by Derickx and Hoeij [18]. Our computation 
suggests that the upper bound b is small, something like 3. 

Let Ki := Q[X]/(P(X)), the Galois representation associated to A(z) is denoted by pi : 
Gal (iff /Q) — > GL^F^), then r(p) mod t is equal to the trace of Frobenius. Since r(p) £ 0(p 6 ), 
take an upper bound L to be O(logp) will be suffice to recover r(p) using Chinese Remainder 
Theorem, the total complexity of getting r(p) is about 

o(^ +26+e ) = 0((logp) 4+Mb+e ) 

£<L, prime 

Remark 4.3. If the constant ui reaches 2.376, the complexity of the algorithm is 0((logp) 8 ' 752+26+e ). 
Moreover, if b is bounded by 3, the complexity is 0((logp) 14 2+e ). 

5. Implementation and results 

The algorithm has been implemented in MAGMA, with a satisfying running time, in view of the 
inefficient operation in the Jacobian, and the calculation of higher level (> 23) becomes feasible. One 
big advantage of the algorithm is that, it is rather straightforward to implement, where the major 
work is dealing with the action of Hecke operators on divisors of function field. The computation 
was done on a computer AMD FX(tm)-6200 Six-Core Processor 3.8GHz. 

For simplicity, we only compute the polynomial corresponding to the projective representation, 
i.e. 

Q(X):= J] (X- 

LeP(Vf) a£L-0 

The simplest case is I = 13, the good prime set is chosen as, 

{19, 23, 29, 43, 53, 61, 67, 71, 79, 83, 89, 109, 127, 149, 157, 163, 179, 193, 211, 223, 229, 233, 239, 241} 

totaly 24 primes, with a product of 52 digits, and 

2535853Q(A) =2535853A 14 + 760835865A 13 + 96570870461 A 12 + 7083218145770A a11 + 

341554192651282A 10 + 11596551892957577X 9 + 288394789072144586A 8 + 
5369247990154339694A 7 + 75509842125272520446A 6 + 
800346109631330635243A 5 + 6303044886777591079517A 4 + 
35793920471135235999031 A 3 + 138667955645963961606844A 2 + 
328650624808255716476451A + 361128579432826593902125 

The computation took several minutes. 

For I = 17, the good prime set consisted of 144 primes with a product of 467 digits, and the 
computation took several hours. 

For i = 19, the good prime needed is about 217 with a product reach 832 digits, and the 
computation took nearly four days. 

For t = 29, for small prime p, the computation of Vt mod p already took nearly one day ( for 
example Vi mod 563 ). The major obstruction is that, operation in the Jacobian is very slow. 

Some verification can be performed, for example using Eichler-Shimura relation to check whether 
the value r(p) mod I is compatible with the factoring type of Q(X) mod p, and check whether 
the discriminant of Q(X) is divided by some large power of £, see[4] for the detail. This all turns 
out to be correct. 
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